function [state algebraic Time]=CICR_CompCore(constants, IC, simulation_duration, Protocol)
global TMAX
TMAX = simulation_duration;
global nAlgebraic;
global shouldCalculateAlgebraic
shouldCalculateAlgebraic = false;
tSpan = [0  TMAX];
options = odeset('RelTol',1e-6);
[Time, state] = ode15s(@(t,y) CICR_func(t,y,constants), tSpan, IC, options);
nAlgebraic = length(Time);
% Algebaraic variabes
shouldCalculateAlgebraic = true;
[junk algebraic] = CICR_func(Time,state',constants);
end